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Abstract. In this work we introduce a new family of twelve-step linear multistep 
methods for the integration of the Schrodinger equation. The new methods are 
constructed by adopting a new methodology which improves the phase lag charac- 
teristics by vanishing both the phase lag function and its first derivatives at a specific 
frequency. This results in decreasing the sensitivity of the integration method on the 
estimated frequency of the problem. The efficiency of the new family of methods is 
proved via error analysis and numerical applications. 
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1. Introduction 

The numerical integration of systems of ordinary differential equa- 
tions with oscillatory solutions has been a subject of research during 
the past decades. This type of ODEs is often met in real problems, 
like the Schrodinger equation and the N-body problem. For problems 
having highly oscillatory solutions, standard methods with unspecial- 
ized use can require a huge number of steps to track the oscillations. 
One way to obtain a more efficient integration process is to construct 
numerical methods with an increased algebraic order, although the 
implementation of high algebraic order meets several difficulties [3] . 

On the other hand, there are some special techniques for optimizing 
numerical methods. Trigonometrical fitting and phase-fitting are some 
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of them, producing methods with variable coefficients, which depend on 
V = ujh, where uj is the dominant frequency of the problem and h is the 
step length of integration. More precisely, the coefficients of a general 
linear method are found from the requirement that it integrates exactly 
powers up to degree p + 1. For problems having oscillatory solutions, 
more efficient methods are obtained when they are exact for every linear 
combination of functions from the reference set 

{l,x,...,x^,e±'^^...,x^e±''-} (1) 

This technique is known as exponential (or trigonometric if /i = ico) 
fitting and has a long history [2], [1]. The set (1) is characterized by 
two integer parameters, K and P . The set in which there is no classical 
component is identified hy K = —1 while the set in which there is no ex- 
ponential fitting component (the classical case) is identified by P = — 1. 
Parameter P will be called the level of tuning. An important property 
of exponential fitted algorithms is that they tend to the corresponding 
classical ones when the involved frequencies tend to zero, a fact which 
allows us to say that exponential fitting represents a natural extension 
of the classical polynomial fitting. The examination of the convergence 
of exponential fitted multistep methods is included in Lyche's theory 
[1] . There is a large number of significant methods presented with high 
practical importance that have been presented in the bibliography. The 
general theory is presented in detail in [5]. 

Considering the accuracy of a method, when solving oscillatory prob- 
lems, it is more appropriate to work with the phase-lag, rather than 
the principal local truncation error. We mention the pioneering paper 
of Brusa and Nigro [4], in which the phase-lag property was intro- 
duced. This is actually another type of a truncation error, i.e. the angle 
between the analytical solution and the numerical solution. On the 
other hand, exponential fitting is accurate only when a good estimate 
of the dominant frequency of the solution is known in advance. This 
means that in practice, if a small change in the dominant frequency is 
introduced, the efficiency of the method can be dramatically altered. 
It is well known that for equations similar to the harmonic oscillator 
the most efficient exponential fitted methods are those with the highest 
tuning level. A lot of significant work has been made during the last 
years in this field, mainly focusing for obvious reasons in the solution 
of the Schrodinger equation (see for example [10]- [132]). 

In this paper we present a new family of methods based on the 12- 
step linear multistep method of Quinlan and Tremaine [6]. The new 
methods are constructed by vanishing the phase-lag function and its 
first derivatives at a predefined frequency. Error analysis and numerical 
experiments show that the new methods exhibit improved characteris- 
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tics concerning the solution of the time-independent Schrodinger equa- 
tion. The paper is organized as follows: In section 2 the general theory 
of the new methodology is presented. In section 3 the new methods 
are described in detail. In section 4 the stability properties of the 
new methods are investigated. Section 5 presents the results from the 
numerical experiments and finally, conclusions are drawn in section 6. 



2. Phase-lag analysis of symmetric multistep methods 

Consider the differential equations 

^ = fit,y), yito) = 2/0, y'{to) = y'o (2) 

and the linear multistep methods 

J J 

Yl y^+j = ^^H ^0 fn+j (3) 

j=0 j=0 

where j/n+i = vik + {n + j)h), fn+j = f{to + in + j)h,y{to + {n + j)h)) 
and h is the step size of the method. With the method (3) , we associate 
the following functional 

J J 
L{h, a, b, y{t)) = ajy{t + j-h)-h^Yl ^oV'^^ + j ' h) (4) 

where o, h are the vectors of coefficients aj and hj respectively, and y{t) 
is an arbitrary function. The algebraic order of the method (3) is p, if 

L(h, a, b, y{t)) = Cj>+2hP+''y^+^\t) + 0{hP+^) (5) 

The coefficients Cq are given 



J 

"'3 



<^0 = Y 

j=0 
J 

j=0 

^ j=0 ' 3=0 

The principal local truncation error (PLTE) is the leading term of (5) 
plte = Cp+2hF+^y^P+'^\t) (7) 
The following assumptions will be considered in the rest of the paper: 
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1. aj = 1, since we can always divide the coefficients of (3) witfi aj. 

2. |ao| + |6o| 7^ 0, since otherwise we can assume that J = J — 1. 

3. J2j=o\^j\ / 0' since otherwise the solution of (3) would be inde- 
pendent of (2). 

4. The method (3) is at least of order one. 

5. The method (3) is zero stable, which means that the roots of the 
polynomial 

J 

p{z) = J2ajz^ (8) 

j=0 

all lie in the unit disc, and those that lie on the unit circle have 
multiplicity one. 

6. The method (3) is symmetric, which means that 

% = bj = bj-j, j = 0(1)J (9) 

It is easily proved that both the order of the method and the step 
number J are even numbers [7]. 

Consider now the test problems 

y"{t) = -cj^yit) (10) 

where a; is a constant. The numerical solution of (10) by applying 
method (3) is described by the difference equation 

J/2 

^ Aj{s'^){yn+j + Vn-j) + Ao{s'^)yn = (11) 

with 

Aj{s^)=ai_. + s^ ■bi_. (12) 
and s = uh. The characteristic equation is then given by 

J/2 

Y,A,{s^){z^ +z-^) + Ao{s') = Q (13) 
i=i 

and the interval of periodicity (0, Sq) is then defined such that for s G 
(0, So) the roots of (13) are of the form 

Zl = e'^^'\ Z2 = e-'^^'\ \zj\ < 1, 3 < j < J (14) 
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where A(s) is a real function of s. The phase-lag PL of the method (3) 
is then defined 

PL = s- X{s) (15) 

and is of order q if 

PL = c- + 0{si+^) (16) 

In general, the coefficients of the method (3) depend on some parameter 
V, thus the coefficients Aj are functions of both and v. The following 
theorem was proved by Simos and Williams [8]: For the symmetric 
method (10) the phase-lag is given 

PL{s, v) = ' (17) 

We are now in position to describe the new methodology. In order to 
efficiently integrate oscillatory problems, it is a good practice to calcu- 
late the coefficients of the numerical method by forcing the phase lag to 
be zero at a specific frequency. But, since the appropriate frequency is 
problem dependent and in general is not always known, we may assume 
that we have an error in the frequency estimation. It would be of great 
importance to force the phase-lag to be insensitive to this error. Thus, 
beyond the vanishing of the phase-lag, we also force its first derivatives 
to be zero. 



3. Construction of the new methods 



3.1. Classical Method 

The family of new methods is based on the 12-step linear multistep 
method of Quinlan and Tremaine [6] which is of the form (3) with 
coefficients 



ao = 1 Ol = —2 02 = 2 03 = —1 

04 = 05 = ae = 

(18) 

h — n h — 90987349 . _ -114798419 t _ 270875723 ^ ^ 

"0 ^ "1 53222400 ""^ 26611200 "3 17740800 

L _ -67855831 l _ 50277247 l _ -253491379 

— 2217600 "5 — 985600 "6 — 44352OO 

The PLTE of the method is given: 

fl630m6103,Mft» 

I 290594304000 \ } ( ^ ^ 
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3.2. New Methods using Phase Fitting 

The methods that are constructed are named as PF-Di, where: 

— PF-DO: the phase lag function is zero at the frequency v = u * h. 

— PF-Dl : the phase lag function and its first derivative are zero at 
the frequency v = uj * h. 

— PF-D2: the phase lag function and its first and second derivatives 

are zero at the frequency v = lo * h. 

— PF-D3: the phase lag function and its first, second and third 
derivatives are zero at the frequency v = uj * h. 

— PF-D4- the phase lag function and its first, second, third and 
fourth derivatives are zero at the frequency v = u * h. 

— PF-D5: the phase lag function and its first, second, third, fourth 
and fifth derivatives are zero at the frequency v = lo * h. 

The coefficients of the methods in the form: 

l.num 



bl 
bl 



bl 



hi 

1, denum 

hi 

2, num 

2, denum 

hi 

3, nuTn 

Ti 

Zjdenum 
hi 

Ti 

i,denum 

hi 
5,nuTn 

Ti 

b,denum 

bi 

o,num 

Ti 

Q,denum 



where the coefficients 6* correspond to the method PF-Di. Since for 
small values ofv, the above formulae are subject to heavy cancelations, 
the Taylor expansions of the coefficients have been calculated as bj^. The 
exact formula of all coefficients are given in appendix. 

The principal local truncation errors of the methods are given by 

16301796103 

pltei = coef fi, 

^ ' 290594304000 ■' ' 
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where 

COe//o= fy(12)^2 ^ y(14)^ ^14 

coefh = Mi0)a;4 + 2y(^^)io ^ + y(^^)) h}^ 

coeff2 = + Sy^^^^u ^ + 3y(i2)^2 + ^(14)^ ^^14 

coe//3 = (y(6)^8 + 4^(8)^6 + g yi^^^co^ + 4y(i2)a;2 + /,14 

coe//4 = (y(^)wio + 5y(6)a; « + 10y(«)a;6 + lOy^^'^'^uo^ 
+5y(i2)^2 + y(14)^ ^14 

coe//5 = (y(2)^i2 + gy(4)^ 10 + i5y(6)^8 + 20y(8)a;6 
+15y(io)a;4 + 6y(i2) ^2 + ^(14)^ ;^14 



4. Stability Analysis 

The stabihty of the new methods is studied by considering the test 
equation 

^ = -a^m (21) 

and the hnear multistep method (3) for the numerical solution. In 
the above equation cr 7^ a; (a; is the frequency at which the phase- 
lag function and its derivatives vanish). By setting s = ah and v = coh, 
we get for the characteristic equation of the applied method 

J/2 

Y,Ms^^)i^'+^~') + Ms',v) = (22) 

where 

Aj{s'^,v) = ai_j{v) + • b^_.{v) (23) 

The motivation of the above analysis is straightforward: Although the 

coefficients of the method (3) are designed in a way that the phase-lag 
and its first derivatives vanish in the frequency cu, the frequency lo itself 
is unknown and only an estimation can be made. Thus, if the correct 
frequency of the problem is a we have to check if the method is stable, 
that is if the roots of the characteristic equation lie in the unit disk. For 
this reason we draw in the s — v plane the areas in which the method 
is stable. Figure 1 shows the stability region for the classical method 
and figure 2 for the six methods (the phase fitted one and those with 
first,second, third, fourth and fifth phase lag derivative elimination). 
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Note here that the s-axis corresponds to the real frequency while the 
?;-axis corresponds to the estimated frequency used to construct the 
parameters of the method. 



5. Numerical Results 

The radial Schrodinger equation is given by: 

/(^) = (^i(i±i) + v{x) - E)y{x)^ (24) 

where ^^^^P is the centrifugal potential, V{x) is the potential, E is the 

Energy and W{x) = ^^^^-P + V{x) is the effective potential. It is valid 
that lirrix^ooy {x) = and therefore limx^ooW {x) = 0. We consider 
that E > and we divide the interval [0, +oo) into subintervals [oj, bi] 
so that W{x) can be considered constant inside each subinterval with 
value Wi- The problem (24) can be expressed now by the equations 

= {Wi - E)y, (25) 

whose solution are 

i{x) = (^yl^e^^^^ + S^e-^^'-^^j (26) 



with Ai,Bi G R. We will integrate problem (24) with ^ = at the 
interval [0, 15] using the well known Woods-Saxon potential: 



Up Uiq ^-^n 



TA/ N "U , >^iy /r,-x 

^(^) = T-r^ + n I „N2 ' Q = e - (27) 



where uq = —50, a = 0.6, Xq = 7, ui = and with boundary 

condition y(0) = 0. The potential V{x) decays more quickly than 
^^^^^ , so for large x (asymptotic region) the Schrodinger equation 
(24) becomes 

y"ix) = (fc^ - E)y{x)^ (28) 

The last equation has two linearly independent solutions kxji{kx) and 
kxni(kx), where ji and n/ arc the spherical Bessel and Neumann func- 
tions and k = \J^^^- — E. When x ^ oo the solution takes the 
asymptotic form 

y{x) ^ Akxji{kx) — Bkxni{kx) 

^ D[sin{kx — vr^) -|- tan{6i) cos{kx — tt^)], (29) 
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where Si is called the scattering phase shift and it is given by the 
following expression: 

= ' (30) 

y{xi+iC{xi) - y{xi)C{xi+i) 

where S{x) = kxji{kx) and C{x) = kxni{kx) and Xi < x^+i and both 
belong to the asymptotic region. Given the energy, wc approximate the 
phase shift, the accurate value of which is ^ for the above problem. 
We will use three different values for the energy: i) 989.701916, ii) 
341.495874 and iii) 163.215341. As for the frequency lo we will use the 
suggestion of Ixaru and Rizea [9] : 



_ j /E^, xG [0,6.5] 
"^'[VE, xe[6.5,15] 

The results are shown in figures 3, 4 and 5. It is clear that the accuracy 
increases as the number of the eliminated derivatives of the phase lag 
function increases. 



6. Conclusions 



We have presented a new family of 12-step symmetric multistep numeri- 
cal methods with improved characteristics concerning the integration of 
the Schrodinger equation. The methods were constructed by adopting 
a new methodology which, except for the phase fitting at a predefined 
frequency, it eliminates the first derivatives of the phase lag function at 
the same frequency. The result is that the phase lag function becomes 
less sensitive on the frequency near the predefined one. This behavior 
compensates the fact that the exact frequency can only be estimated. 
Experimental results demonstrate this behavior by showing that the ac- 
curacy is increased as the number of the derivatives that are eliminated 
is increased. 



References 

1. T. Lyche, Chebyshevian multistep methods for Ordinary DifTerential Eqa- 
tions, Num. Math. 19, 65-75 (1972) 

2. Gautschi, W.: Numerical integration of ordinary differential equations based 
on trigonometric polynomials. Numer. Math. 3, 381397 (1961) 

3. Quinlan, G.: Resonances and instabilities in symmetric multistep methods, 
preprint arXiv astro-ph/9901136 (1999) 

4. L. Brusa, L.N.: A one-step method for direct integration of structural dynamic 
equations. Int. J. Num. Methods Engrg. 15, 685699 (1980) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 9 



10 D.S. Vlachos, Z.A. Anastassi, T.E. Simos 

5. L.Gr. Ixaru, G.V.Berghe: Kluwer Academic Publishers, Dor- 
drecht /Boston/London (2004) 

6. Quinlan, D., Tremaine, S.: Symmetric multistep methods for the numerical 
integration of planetary orbits. The Astronomical Journal 100(5), 16941700 
(1990) 

7. Lambert, J., Watson, L: Symmetric multistep methods for periodic initial 
values problems. J. Inst. Math. Appl. 18, 189202 (1976) 

8. Simos, T., Williams, P.: On finite difference methods for the solution of the 
Schrodinger equation. Comput. Chem. 23, 513554 (1999) 

9. Ixaru, L.G., Rizea, M.: Comparison of some four-step methods for the numer- 
ical solution of the Schrodinger equation. Computer Physics Communications 
38(3), 329337 (1985) 

10. L.Gr. Ixaru and M. Micu, Topics in Theoretical Physics. Central Institute of 
Physics, Bucharest, 1978. 

11. L.D. Landau and F.M. Lifshitz: Quantum Mechanics. Pergamon, New York, 
1965. 

12. I. Prigogine, Stuart Rice (Eds): Advances in Chemical Physics Vol. 93: New 
Methods in Computational Quantum Mechanics, John Wiley & Sons, 1997. 

13. G. Herzberg, Spectra of Diatomic Molecules, Van Nostrand, Toronto, 1950. 

14. T.E. Simos, Atomic Structure Computations in Chemical Modelling: Appli- 
cations and Theory (Editor: A. Hinchliffe, UMIST), The Royal Society of 
Chemistry 38-142(2000). 

15. T.E. Simos, Numerical methods for ID, 2D and 3D differential equations 
arising in chemical problems. Chemical Modelling: Application and Theory, 
The Royal Society of Chemistry, 2(2002), 170-270. 

16. T.E. Simos: Numerical Solution of Ordinary Differential Equations with Pe- 
riodical Solution. Doctoral Dissertation, National Technical University of 
Athens, Greece, 1990 (in Greek). 

17. A. Konguetsof and T.E. Simos, On the Construction of exponentially-fitted 
methods for the numerical solution of the Schrodinger Equation, Journal of 
Computational Methods in Sciences and Engineering 1 143-165(2001). 

18. A.D. Raptis and A.C. Allison: Exponential - fitting methods for the numerical 
solution of the Schrodinger equation. Computer Physics Communications, 14 
1-5(1978). 

19. A.D. Raptis, Exponential multistep methods for ordinary differential equa- 
tions. Bull. Greek Math. Soc. 25 113-126(1984). 

20. L.Gr. Ixaru, Numerical Methods for Differential Equations and Applications, 
Reidel, Dordrecht - Boston - Lancaster, 1984. 

21. L.Gr. Ixaru and M. Rizea, A Numerov-like scheme for the numerical solution 
of the Schrodinger equation in the deep continuum spectrum of energies. 
Comput. Phys. Commun. 19 23-27(1980). 

22. T. E. Simos, P. S. Williams: A New Runge-Kutta-Nystrom Method with 
Phase-Lag of Order Infinity for the Numerical Solution of the Schrodinger 
Equation, MATCH Commun. Math. Comput. Chem. 45 123-137(2002). 

23. T. E. Simos, Multidorivativc Methods for the Numerical Solution of the 
Schrodinger Equation, MATCH Commun. Math. Comput. Chem. 45 7- 
26(2004). 

24. A.D. Raptis, Exponentially-fitted solutions of the eigenvalue Shrodingcr equa- 
tion with automatic error control. Computer Physics Communications, 28 
427-431(1983) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 10 



High Order Multistep Methods with Improved Phase-Lag Characteristics 11 

25. A.D. Raptis, On the numerical solution of the Schrodinger equation, Computer 
Physics Communications, 24 1-4(1981) 

26. Zacharoula Kalogiratou and T.E. Simos, A P-stable exponentially-fitted 
method for the numerical integration of the Schrodinger equation, Applied 
Mathematics and Computation, 112 99-112(2000). 

27. A.D. Raptis and T.E. Simos, A four-step phase-fitted method for the numer- 
ical integration of second order initial-value problem, BIT, 31 160-168(1991). 

28. Peter Henrici, Discrete variable methods in ordinary differential equations, 
John Wiley & Sons, 1962. 

29. M.M. Chawla, Uncoditionally stable Noumerov-type methods for second order 
differential equations, BIT, 23 541-542(1983). 

30. M. M. Chawla and P. S. Rao, A Noumorov-typc method with minimal phase- 
lag for the integration of second order periodic initial- value problems. Journal 
of Computational and Applied Mathematics 11(3) 277-281(1984) 

31. Z.A. Anastassi, T.E. Simos, A family of cxponcntially-fittcd Rungc-Kutta 
methods with exponential order up to three for the numerical solution of the 
Schrodinger equation, J. Math. Chem 41 (1) 79-100 (2007) 

32. T. Monovasilis, Z. Kalogiratou , T.E. Simos, Trigonometrically fitted and 
exponentially fitted symplectic methods for the numerical integration of the 
Schrodinger equation, J. Math. Chem 40 (3) 257-267 (2006) 

33. G. Psilioyios, T.E. Simos, The numerical solution of the radial Schrodinger 
equation via a trigonometrically fitted family of seventh algebraic order 
Predictor-Corrector methods, J. Math. Chem 40 (3) 269-293 (2006) 

34. T.E. Simos, A four-step exponentially fitted method for the numerical solution 
of the Schrodinger equation, J. Math. Chem 40 (3) 305-318 

(2006) 

35. T. Monovasilis, Z. Kalogiratou , T.E. Simos, Exponentially fitted symplectic 
methods for the numerical integration of the Schrodinger equation J. Math. 
Chem 37 (3) 263-270 (2005) 

36. Z. Kalogiratou , T. Monovasilis, T.E. Simos, Numerical solution of the 
two-dimensional time independent Schrodinger equation with Numerov-type 
methods J. Math. Chem 37 (3) 271-279 (2005) 

37. Z.A. Anastassi, T.E. Simos, Trigonometrically fitted Runge-Kutta methods 
for the numerical solution of the Schrodinger equation J. Math. Chem 37 (3) 
281-293 (2005) 

38. G. Psihoyios, T.E. Simos, Sixth algebraic order trigonometrically fit- 
ted predictor-corrector methods for the numerical solution of the radial 
Schrodinger 

equation, J. Math. Chem 37 (3) 295-316 (2005) 

39. D.P. Sakas, T.E. Simos, A family of multiderivative methods for the numerical 
solution of the Schrodinger equation, J. Math. Chem 37 (3) 

317-331 (2005) 

40. T.E. Simos, Exponentially - fitted multiderivative methods for the numerical 
solution of the Schrodinger equation, J. Math. Chem 36 (1) 13-27 (2004) 

41. K. Tsclios, T.E. Simos, Symplectic methods of fifth order for the numerical 
solution of the radial Shrodinger equation, J. Math. Chem 35 (1) 

55-63 (2004) 

42. T.E. Simos, A family of trigonomctrically-fittcd symmetric methods for the 
efficient solution of the Schrodinger equation and related problems J. Math. 
Chem 34 (1-2) 39-58 JUL 2003 



QT12_JMC.tex; 15/11/2008; 10:30; p. 11 



12 D.S. Vlachos, Z.A. Anastassi, T.E. Simos 

43. K. Tselios, T.E. Simos, Symplectic methods for the numerical solution of the 
radial Shrodingor equation, J, Math. Chem 34 (1-2) 83-94 (2003) 

44. J. Vigo-Aguiar J, T.E. Simos, Family of twelve steps exponential fitting 
symmetric multistep methods for the numerical solution of the Schrodinger 
equation, J. Math. Chem 32 (3) 257-270 (2002) 

45. G. Avdelas, E. Kefalidis, T.E. Simos, New P-stable eighth algebraic order 
exponontially-fittod methods for the numerical integration of the Schrodinger 
equation, J. Math. Chem 31 (4) 371-404 (2002) 

46. T.E. Simos, J. Vigo-Aguiar, Symmetric eighth algebraic order methods with 
minimal phase-lag for the numerical solution of the Schrodinger equation J. 
Math. Chem 31 (2) 135-144 (2002) 

47. Z. Kalogiratou , T.E. Simos, Construction of trigonometrically and exponen- 
tially fitted Runge-Kutta-Nystrom methods for the numerical solution of the 
Schrodinger equation and related problems a method of 8th algebraic order, 
J. Math. Chem 31 (2) 211-232 

48. T.E. Simos, J. Vigo-Aguiar, A modified phase-fitted Runge-Kutta method for 
the numerical solution of the Schrodinger equation, J. Math. Chem 30 (1) 
121-131 (2001) 

49. G. Avdelas, A. Konguetsof, T.E. Simos, A generator and an optimized gen- 
erator of high-order hybrid explicit methods for the numerical solution of the 
Schrodinger equation. Part 1. Development of the basic method, J. Math. 
Chem 29 (4) 281-291 (2001) 

50. G. Avdelas, A. Konguetsof, T.E. Simos, A generator and an optimized gen- 
erator of high-order hybrid explicit methods for the numerical solution of the 
Schrodinger equation. Part 2. Development of the generator; optimization of 
the generator and numerical results, J. Math. Chem 29 (4) 293-305 (2001) 

51. J. Vigo-Aguiar, T.E. Simos, A family of P-stable eighth algebraic order 
methods with exponential fitting facilities, J. Math. Chem 29 (3) 177-189 
(2001) 

52. T.E. Simos, A new explicit Bessel and Neumann fitted eighth algebraic order 
method for the numerical solution of the Schrodinger equation J. Math. Chem 
27 (4) 343-356 (2000) 

53. G. Avdelas, T.E. Simos, Embedded eighth order methods for the numerical 
solution of the Schrodinger equation, J. Math. Chem 26 (4) 327-341 1999, 

54. T.E. Simos, A family of P-stable exponentially-fitted methods for the nu- 
merical solution of the Schrodinger equation, J. Math. Chem 25 (1) 65-84 
(1999) 

55. T.E. Simos, Some embedded modified Runge-Kutta methods for the numerical 
solution of some specific Schrodinger equations, J. Math. Chem24 (1-3) 23-37 

(1998) 

56. T.E. Simos, Eighth order methods with minimal phase-lag for accurate com- 
putations for the elastic scattering phase-shift problem, J. Math. Chem 21 
(4) 359-372 (1997) 

57. P. Amodio, I. Gladwell and G. Romanazzi, Numerical Solution of General 
Bordered ABD Linear Systems by Cyclic Reduction, JNAIAM J. Numer. 
Anal. Indust. Appl. Math 1(1) 5-12(2006) 

58. S. D. Capper, J. R. Cash and D. R. Moore, Lobatto-ObrechkofF Formulae for 
2nd Order Two-Point Boundary Value Problems, JNAIAM J. Numer. Anal. 
Indust. Appl. Math 1(1) 13-25 (2006) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 12 



High Order Multistep Methods with Improved Phase-Lag Characteristics 13 

59. S. D. Capper and D. R. Moore, On High Order MIRK Schemes and Hermite- 
Birkhoff Interpolants, JNAIAM J. Numer. Anal. Indust. Appl. Math 1(1) 
27-47 (2006) 

60. J. R. Cash, N. Sumaxti, T. J. Abdulla and I. Vieira, The Derivation of In- 
terpolants for Nonlinear Two-Point Boundary Value Problems, JNAIAM J. 
Numer. Anal. Indust. Appl. Math 1(1) 49-58 (2006) 

61. J. R. Cash and S. Girdlcstonc, Variable Step Rungo-Kutta-Nystrm Methods 
for the Numerical Solution of Reversible Systems, JNAIAM J. Numer. Anal. 
Indust. Appl. Math 1(1) 59-80 (2006) 

62. Jeff R. Cash and Prancesca Mazzia, Hybrid Mesh Selection Algorithms 
Based on Conditioning for Two-Point Boundary Value Problems, JNAIAM 
J. Numer. Anal. Indust. Appl. Math 1(1) 81-90 (2006) 

63. Felice lavernaro, Francesca Mazzia and Donato Trigiante, Stability and Condi- 
tioning in Numerical Analysis, JNAIAM J. Numer. Anal. Indust. Appl. Math 
1(1) 91-112 (2006) 

64. Felice lavernaro and Donato Trigiante, Discrete Conservative Vector Fields 
Induced by the Trapezoidal Method, JNAIAM J. Numer. Anal. Indust. Appl. 
Math 1(1) 113-130 (2006) 

65. Prancesca Mazzia, Alessandra Sestini and Donato Trigiante, BS Linear Mul- 
tistep Methods on Non-uniform Meshes, JNAIAM J. Numer. Anal. Indust. 
Appl. Math 1(1) 131-144 (2006) 

66. L.F. Shampine, P.H. Muir, H. Xu, A User-Friendly Fortran BVP Solver, 
JNAIAM J. Numer. Anal. Indust. Appl. Math 1(2) 201-217 (2006) 

67. G. Vanden Berghe and M. Van Daele, Exponentially- fitted Strmer/Verlet 
methods, JNAIAM J. Numer. Anal. Indust. Appl. Math 1(3) 241-255 (2006) 

68. L. Aceto, R. Pandolfi, D. Trigiante, Stability Analysis of Linear Multistep 
Methods via Polynomial Type Variation, JNAIAM J. Numer. Anal. Indust. 
Appl. Math 2(1-2) 1-9 (2007) 

69. G. Psihoyios, A Block Implicit Advanced Step-point (BIAS) Algorithm for 
Stiff Differential Systems, Computing Letters 2(1-2) 51-58(2006) 

70. W.H. Enright, On the use of 'arc length' and 'defect' for mesh selection for 
differential equations. Computing Letters 1(2) 47-52(2005) 

71. T.E. Simos, P-stable Four-Step Exponentially-Fitted Method for the Nu- 
merical Integration of the Schrodinger Equation, Computing Letter 1(1) 
37-45(2005). 

72. T.E. Simos, Stabilization of a Four-Step Exponentially-Fitted Method and 
its Application to the Schrodinger Equation, International Journal of Modem 
Physics C 18(3) 315-328(2007). 

73. Zhongcheng Wang, P-stable linear symmetric multistep methods for periodic 
initial-value problems, Computer Physics Communications 171 162174(2005) 

74. T.E. Simos, A Runge-Kutta Fehlberg method with phase-lag of order in- 
finity for initial value problems with oscillating solution. Computers and 
Mathematics with Applications 25 95-101(1993). 

75. T.E. Simos, Runge-Kutta interpolants with minimal phase-lag. Computers 
and Mathematics with Applications 26 43-49(1993). 

76. T.E. Simos, Runge-Kutta-Nystrom interpolants for the numerical integra- 
tion of special second-order periodic initial-value problems. Computers and 
Mathematics with Apphcoiions 26 7-15(1993). 

77. T.E. Simos and G.V. Mitsou, A family of four-step exponential fitted methods 
for the numerical integration of the radial Schrodinger equation, Computers 
and Mathematics with Applications 28 41-50(1994). 



QT12_JMC.tex; 15/11/2008; 10:30; p. 13 



14 D.S. Vlachos, Z.A. Anastassi, T.E. Simos 

78. T.E. Simos and G. Mousaxlis, A two-step method for the numerical solution of 
the radial Schrdinger equation, Computers and Mathematics with Applications 
29 31-37(1995). 

79. G. Avdelas and T.E. Simos, Block Runge-Kutta methods for periodic 
initial-value problems. Computers and Mathematics with Applications 31 69- 
83(1996). 

80. G. Avdelas and T.E. Simos, Embedded methods for the numerical solution of 
the Schrodinger equation. Computers and Mathematics with Applications 31 
85-102(1996). 

81. G. Papakaliatakis and T.E. Simos, A new method for the numerical solution 
of fourth order BVPs with oscillating solutions, Computers and Mathematics 
with Applications 32 1-6(1996). 

82. T.E. Simos, An extended Numerov-type method for the numerical solution of 
the Schrodinger equation. Computers and Mathematics with Applications 33 
67-78(1997). 

83. T.E. Simos, A new hybrid imbedded variable-step procedure for the numerical 
integration of the Schrodinger equation. Computers and Mathematics with 
Applications 36 51-63(1998). 

84. T.E. Simos, Bessel and Neumann Fitted Methods for the Numerical Solution 
of the Schrodinger equation, Computers & Mathematics with Applications 42 
833-847(2001). 

85. A. Konguetsof aud T.E. Simos, An exponentially-fitted and trigonometrically- 
fitted method for the numerical solution of periodic initial-value problems. 
Computers and Mathematics with Applications 45 547-554(2003). 

86. Z.A. Anastassi and T.E. Simos, An optimized Runge-Kutta method for 
the solution of orbital problems, Journal of Computational and Applied 
Mathematics 175(1) 1-9(2005) 

87. G. Psihoyios and T.E. Simos, A fourth algebraic order trigonometrically fit- 
ted predictor-corrector scheme for IVPs with oscillating solutions, Journal of 
Computational and Applied Mathematics 175(1) 137-147(2005) 

88. D.P. Sakas and T.E. Simos, Multiderivative methods of eighth algrebraic order 
with minimal phase-lag for the numerical solution of the radial Schrodinger 
equation. Journal of Computational and Applied Mathematics 175(1) 161- 
172(2005) 

89. K. Tselios and T.E. Simos, Runge-Kutta methods with minimal dispersion 
and dissipation for problems arising from computational acoustics. Journal of 
Computational and Applied Mathematics 175(1) 173-181(2005) 

90. Z. Kalogiratou and T.E. Simos, Newton-Cotes formulae for long-time integra- 
tion. Journal of Computational and Applied Mathematics 158(1) 75-82(2003) 

91. Z. Kalogiratou, T. Monovasilis and T.E. Simos, Symplectic integrators for the 
numerical solution of the Schrodinger equation. Journal of Computational and 
Applied Mathematics 158(1) 83-92(2003) 

92. A. Konguetsof and T.E. Simos, A generator of hybrid symmetric four-step 
methods for the numerical solution of the Schrodinger equation. Journal of 
Computational and Applied Mathematics 158(1) 93-106(2003) 

93. G. Psihoyios and T.E. Simos, Trigonometrically fitted predictor-corrector 
methods for IVPs with oscillating solutions. Journal of Computational and 
Applied Mathematics 158(1) 135-144(2003) 

94. Ch. Tsitouras and T.E. Simos, Optimized Runge-Kutta pairs for problems 
with oscillating solutions. Journal of Computational and Applied Mathematics 
147(2) 397-409(2002) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 14 



High Order Multistep Methods with Improved Phase-Lag Characteristics 15 

95. T.E. Simos, An exponentially fitted eighth-order method for the numerical 
solution of the Schrodinger equation, Journal of Computational and Applied 
Mathematics 108(1-2) 177-194(1999) 

96. T.E. Simos, An accurate finite difference method for the numerical solution of 
the Schrodinger equation. Journal of Computational and Applied Mathematics 
91(1) 47-61(1998) 

97. R.M. Thomas and T.E. Simos, A family of hybrid exponentially fitted 
predictor-corrector methods for the numerical integration of the radial 
Schrodinger equation. Journal of Computational and Applied Mathematics 
87(2) 215-226(1997) 

98. Z.A. Anastassi and T.E. Simos: Special Optimized Runge-Kutta methods for 
IVPs with Oscillating Solutions, International Journal of Modern Physics C, 
15, 1-15 (2004) 

99. Z.A. Anastassi and T.E. Simos: A Dispersive-Fitted and Dissipative-Fitted 
Explicit Rungc-Kutta method for the Numerical Solution of Orbital Problems, 
New Astronomy, 10, 31-37 (2004) 

100. Z.A. Anastassi and T.E. Simos: A Trigonomctrically-Fitted Runge-Kutta 
Method for the Numerical Solution of Orbital Problems, New Astronomy, 
10, 301-309 (2005) 

101. T.V. TriantafyUidis, Z.A. Anastassi and T.E. Simos: Two Optimized Runge- 
Kutta Methods for the Solution of the Schr?dinger Equation, MATCH 
Comnmn. Math. Comput. Chem., 60, 3 (2008) 

102. Z.A. Anastassi and T.E. Simos, Trigonometrically Fitted Fifth Order Runge- 
Kutta Methods for the Numerical Solution of the Schrodinger Equation, 
Mathematical and Computer Modelling, 42 (7-8), 877-886 (2005) 

103. Z.A. Anastassi and T.E. Simos: New Trigonometrically Fitted Six-Step Sym- 
metric Methods for the Efficient Solution of the Schrodinger Equation, 
MATCH Commun. Math. Comput. Chem., 60, 3 (2008) 

104. G.A. Panopoulos, Z.A. Anastassi and T.E. Simos: Two New Optimized 
Eight-Step Symmetric Methods for the Efficient Solution of the Schrodinger 
Equation and Related Problems, MATCH Commun. Math. Comput. Chem., 
60, 3 (2008) 

105. Z.A. Anastassi and T.E. Simos: A Six-Step P-stable Trigonometrically-Fitted 
Method for the Numerical Integration of the Radial Schrodinger Equation, 
MATCH Commun. Math. Comput. Chem., 60, 3 (2008) 

106. Z.A. Anastassi and T.E. Simos, A family of two-stage two-step methods for 
the numerical integration of the Schrodinger equation and related IVPs with 
oscillating solution. Journal of Mathematical Chemistry, Article in Press, 
Corrected Proof 

107. T.E. Simos and P.S. Williams, A finite-difference method for the numerical 
solution of the Schrodinger equation. Journal of Computational and Applied 
Mathematics 79(2) 189-205(1997) 

108. G. Avdelas and T.E. Simos, A generator of high-order embedded P-stable 
methods for the numerical solution of the Schrodinger equation. Journal of 
Computational and Applied Mathematics 72(2) 345-358(1996) 

109. R.M. Thomas, T.E. Simos and G.V. Mitsou, A family of Numerov-type ex- 
ponentially fitted predictor-corrector methods for the numerical integration 
of the radial Schrodinger equation. Journal of Computational and Applied 
Mathematics 67(2) 255-270(1996) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 15 



16 D.S. Vlachos, Z.A. Anastassi, T.E. Simos 

110. T.E. Simos, A Family of 4-Step Exponentially Fitted Predictor-Corrector 
Methods for the Numerical-Integration of The Sclirodingor-Equation, Journal 
of Computational and Applied Mathematics 58(3) 337-344(1995) 

111. T.E. Simos, An Explicit 4-Step Phase-Fitted Method for the Numerical- 
Integration of 2nd-order Initial- Value Problems, Journal of Computational 
and Applied Mathematics 55(2) 125-133(1994) 

112. T.E. Simos, E. Dimas and A.B. Sideridis, A Runge-Kutta-Nystrom Method 
for the Numerical-Integration of Special 2nd-order Periodic Initial- Value 
Problems, Journal of Computational and Applied Mathematics 51(3) 317- 
326(1994) 

113. A.B. Sideridis and T.E. Simos, A Low-Order Embedded Runge-Kutta Method 
for Periodic Initial- Value Problems, Journal of Computational and Applied 
Mathematics 44(2) 235-244(1992) 

114. T.E. Simos amd A.D. Raptis, A 4th-order Bessel Fitting Method for the 
Numerical-Solution of the SchrOdinger-Equation, Journal of Computational 
and Applied Mathematics 43(3) 313-322(1992) 

115. T.E. Simos, Explicit 2-Step Methods with Minimal Phase-Lag for the 
Numerical-Integration of Special 2nd-order Initial- Value Problems and their 
Application to the One-Dimensional Schrodinger-Equation, Journal of Com- 
putational and Applied Mathematics 39(1) 89-94(1992) 

116. T.E. Simos, A 4-Step Method for the Numerical-Solution of the Schrodinger- 
Equation, Journal of Computational and Applied Mathematics 30(3) 251- 
255(1990) 

117. CD. Papageorgiou, A.D. Raptis and T.E. Simos, A Method for Com- 
puting Phase-Shifts for Scattering, Journal of Computational and Applied 
Mathematics 29(1) 61-67(1990) 

118. A.D. Raptis, Two-Step Methods for the Numerical Solution of the Schrodinger 
Equation, Computing 28 373-378(1982). 

119. T.E. Simos. A new Numerov-type method for computing eigenvalues and res- 
onances of the radial Schrodinger equation. International Journal of Modern 
Physics C-Physics and Computers, 7(1) 33-41(1996) 

120. T.E. Simos, Predictor Corrector Phase-Fitted Methods for Y" =F(X,Y) and an 
Application to the Schrodinger-Equation, International Journal of Quantum 
Chemistry, 53(5) 473-483(1995) 

121. T.E. Simos, Two-step almost P-stable complete in phase methods for the 
numerical integration of second order periodic initial-value problems. Inter. 
J. Comput. Math. 46 77-85(1992). 

122. R. M. Corless, A. Shakoori, D.A. Aruliali, L. Gonzalez- Vega, Barycentric 
Hermite Interpolants for Event Location in Initial- Value Problems, JNAIAM 
J. Nurner. Anal. Indust. Appl. Math, 3, 1-16 (2008) 

123. M. Dewar, Embedding a General-Purpose Numerical Library in an Interactive 
Environment, JNAIAM J. Numer. Anal. Indust. Appl. Math, 3, 17-26 (2008) 

124. J. Kierzenka and L.F. Shampine, A BVP Solver that Controls Residual and 
Error, JNAIAM J. Numer. Anal. Indust. Appl. Math, 3, 27-41 (2008) 

125. R. Knapp, A Method of Lines Framework in Mathematica, JNAIAM J. 
Numer. Anal. Indust. Appl. Math, 3, 43-59 (2008) 

126. N. S. Nedialkov and J. D. Pryce, Solving Differential Algebraic Equations 
by Taylor Series (III): the DAETS Code, JNAIAM J. Numer. Anal. Indust. 
Appl. Math, 3, 61-80 (2008) 

127. R. L. Lipsman, J. E. Osborn, and J. M. Rosenberg, The SCHOL Project at 
the University of Maryland: Using Mathematical Software in the Teaching of 



QT12_JMC.tex; 15/11/2008; 10:30; p. 16 



High Order Multistep Methods with Improved Phase-Lag Characteristics 17 

Sophomore Differential Equations, JNAIAM J. Numer. Anal. Indust. Appl. 
Math, 3, 81-103 (2008) 

128. M. Sofroniou and G. Spaletta, Extrapolation Methods in Mathematica, 
JNAIAM J. Numer. Anal. Indust. Appl. Math, 3, 105-121 (2008) 

129. R. J. Spiteri and Thian-Peng Ter, pythNon: A PSE for the Numerical Solution 
of Nonlinear Algebraic Equations, JNAIAM J. Numer. Anal. Indust. Appl. 
Math, 3, 123-137 (2008) 

130. S.P. Corwin, S. Thompson and S.M. White, Solving ODEs and DDEs with 
Impulses, JNAIAM J. Numer. Anal. Indust. Appl. Math, 3, 139-149 (2008) 

131. W. Weckesser, VFGEN: A Code Generation Tool, JNAIAM J. Numer. Anal. 
Indust. Appl. Math, 3, 151-165 (2008) 

132. A. Wittkopf, Automatic Code Generation and Optimization in Maple, 
JNAIAM J. Numer. Anal. Indust. Appl. Math, 3, 167-180 (2008) 



QT12_JMC.tex; 15/11/2008; 10:30; p. 17 



18 D.S. Vlachos, Z.A. Anastassi, T.E. Simos 

Appendix 



Method PF-DO: 



h\^^^^„^ = ((-124184636 cos (t;) + 70378348 cos (2i;) 
-24862148 cos(3v) + 5153611 cos(4v))t>2 
+25(3013169i;2 - 16128 cos(3u) 

+32256 cos (4i;) - 32256 cos(5t;) + 16128 cos(6z;))) csc^°(|) 
^2,n«m = (1^^(159588050008(1;) - 85350160 cos (2^;) 
+16708985 cos (3t;) - 5153611 cos (5t;)) - 16(6496079t;2 
-252000 cos (3w) + 504000 cos(4z;) 
-504000 cos(5f ) + 252000 cos(6z;))) csci°(|) 
bl,num = ((-367257540 COS (t;) + 183567900 cos (2?;) 
-16708985 cos (4u) + 24862148 cos (5w))?;2 
+81(3175117t;2 - 224000 cos (3v) + 448000 cos (4t;) 
-448000 cos (5t;) + 224000 cos (6?;))) esc i°(f) 
i>Xnum = (30675810 cos (i;)i;2 - 429587881;^ 
+75(161280 - 611893^2) cos(3w) 
+140(152411t;2 - 172800) cos (4z;) 

+(24192000 - 17594587^2) cos(5i;) - 12096000 cos(6i;)) csc^O(|) 

= (-61351620 cos (2z;)z;2 + 85936557i;2 
+30(6120959t;2 - 1411200) cos(3t;) 
+25(3386880 - 3191761t;2) cos(4t;) 

+(62092318?;2 - 84672000) cos(5i;) + 42336000 cos (6?;)) csc^°(|) 
bl,num = ((-171873114cos(?;) + 171835152 cos (2i;) 
-257184477 cos (3v) + 103937264 cos (4t;) 
-75329225 cos(5i;))?;2 + 50803200(cos(3t;) - 2cos(4i;) 
+2cos(5u) - cos(6i;)))csci°(|) 
b%enom = 206438400t;2 
^Idenom = 206438400^;^ 
bldenom = 206438400z;2 
bldenam = 516096001;^ 
bUnom = 103219200t;2 
bldenom = 103219200x;2 



(32) 
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_ 90987349 16301796103 v'^ 2012122579?;^ 

'^'^ ~ 53222400 ~ 290594304000 ^ 581188608000 
48410309 (''^ 991915331743 
+ 



1184440320000 2838385676206080000 
__ 114798419 16301796103 _ 2012122579?;'' 

~ ~ 26611200 ^ 29059430400 „ " 58118860800 
48410309 991945331743 

+ ... 



118444032000 283838567620608000 
^0 _ 270875723 16301796103 r- 2012122579?/^ 

'^'^ ~ 17740800 6457651200 , ^ 12915302400 
48410309 991945331743 

"26320896000 ^ 63075237249024000 ~ ' ' ' 

^0 __ 67855831 16301796103 ?;^ 2012122579?;^ (33) 

^'^ ~ ~ 2217600 ^ 2421619200_ 4843238400 
48410309 991945331713 

+ ... 



9870336000 23653213968384000 
_ 50277247 16301796103 2012122579?;^ 

'^'^ ~ 985600 1383782400 „ 2767564800 

48410309 991945331743 

~ 5640192000 ^ 13516122267648000 ~ 
_ 253491379 16301790103 2012122579?;'^ 

'^'^ ~ 4435200 ^ 1153152000 2306304000 
48410309 991945331743 

4700160000 ~ 11263435223040000 " ' 
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Method PF-Dl: 

Knnm = csci°(|) sec(| ) (53760(2 COS (v) + 2cos(2i;) + 2cos(3?;) 

+2cos(5 v) + l)sin''^(|) + z;((304823cos (|) - 828603 cos(|^) 
+554639 cos(f) - 272779 cos C-^))v'^ + 6720(7cos(^) 
-15cos(f ) + 2(9cos( f ) - 6cos(ii^^) +COS {^))))) 
barium = csc^d )(w((393310 cos(2w) - 211395 COS (3?;) 
+272779 cos(4i;))7;2 ^ (122851 v'^ + 1088640) cos(t;) 
-15(-1817f2 + 74368 cos(2u) - 65856 cos(3i;) + 45248 cos(4v) 
-19040 cos (5z;) + 1344 cos (6w) + 224 cos (7t;) + 36288)) csc^d) 
-26880(10 cos(t;) + 10cos(2t;) + 9cos(3i;) + 2cos(4i;) 
+8cos(5z;) +cos(6 t;) + 5) sec(|)) 

H,num = csc^(^)(26880(90cos (v) + 90cos(2t;) + 74cos(3i;) 

+32cos(4t;) + 58cos(5i;) + 16 cos {6v) + 45)sec(|) 

— y((2519110cos(2 v) - 654530 cos (3t;) + 272779(6 cos (4t>) 

+ cos(5?;))) + (1288471^2 + 4717440) cos(^;) + 8(47587 

-643440 cos(2i;) + 593880 cos (3i;) - 416640 cos (4t;) 

+159600 cos (5u) + 18480 cos(6i;) 

-6720 cos(7i;) - 294840)) CSC 3(|)) 

bl,num = csc^(|)(t;(1578002 v"^ + 2(161489^;^ + 725760) cos(t;) 
+160(20521 - 13776) cos(2t>) + 5(478464 - m9v^) cos(3 v) 
+2(738791f 2 - 900480) cos(4^;) + 3(205341 v"^ + 156800) cos(5?;) 
+510720 cos(6w) - 94080 cos(7 v) - 725760) csc^d) 
-107520(30 cos (v) + 30cos(2t;) + 23cos(3i;) + 14cos(4i;) 
+16cos(5t;) + 7cos(6 v) + 15)sec(f )) 

bl^^^^ = - csc^d )(7;(923342 v'^ + 28(3433391;^ - 181440) cos(t;) 

+20(240593 + 103488) cos(2t>) + 5(832859^;^ + 18816) cos(3 v) 

+2(1574159?;^ - 376320) cos(4t;) + 3(913431 v'^ - 313600) cos(5f ) 

+2446080 cos (6z;) - 376320 cos(7 v) + 2540160) csc^d) 

-376320(30 cos (v) + 30cos(2?;) + 22cos(3i;) + 16cos(4i;) 

+14 cos(5?;) + 8 cos(6 v) + 15) sec(|)) 

blnum = csc^d )(v(7 (65953?;^ + 181440) + 7(379709^;^ 

-362880) cos {v) + 2(694753^2 + 799680) cos(2t>) + 4(360653 v'^ 

-199920) cos(3f ) + 15(52889t;2 + 21952) cos(4 v) 

+(875393i;2 - 540960) cos (5?;) + 799680 cos (6?;) 

-117600 cos(7?;))csc3(|) - 188160(18 cos (v) + 18cos(2?;) 

+13 cos(3t;) + 10 cos(4t;) + 8 cos(5i;) + 5 cos(6 v) + 9) sec(|)) 

bUenom = 6881280^;^ 

bidenom = 1720320^;=^ 

bUnom = MmAOv' 
bidenom = ™320v^ 



b\,denom = 3440640^;^ 

bldenom = 860160^3 



(34) 
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^1 _ 90987349 16301796103 1532031563 t;^ 

^'^ ~ 53222400 T 145297152000 ^ 268240896000 
31987133939 5466168990203 

+ 



592812380160000 2838385676206080000 
1 _ 114798419 16301796103 28198975459 

'^'^ ~ ~ 26611200 ^ 14529715200 249080832000 
88492028011 533071354889581t;^ 

"^11856247603200 ~ 1419192838103040000 " ' 
1 _ 270875723 _ 16301796103 117200518583 t;^ 

^■^ ~ 17740800 , 3228825600 ^ 166053888000 
_ 2285150767769 308808361613933 

39520825344000 ^ 105125395415040000 ~ ' ' ' 
^1 __ 67855831 16301796103 t;2 36423021893 t;^ 

^'^ ~ ~ 2217600 1210809600 16144128000 
989757481543 1207392034807669 



4940103168000 118266069841920000 
1 _ 50277247 16301796103 360410789243 

'^'^ ~ 985600 , 691891200 ^ 83026944000 
25035151487 125088382253381 ?/ 2919882783787129 v^^ 

62731468800 6143691939840000 ~ 4460320348323840000 
1 _ 253491379 16301796103 222767191987 

^'^ ~ 4435200 576576000 „ 41513472000 

3516569786117 95477706318691 

^ 7057290240000 3754478407680000 ^ ' ' ' 

(35) 
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Method PF-D2: 

Knum = -30720(i;2 - 1) cos{v) - 10(-387csc ^(|) 

+72(128 cos(7;) + 127) csc'^{v) + 2048 (v^ - 3)) - {8{-47v'^ 
-2040v^ - (40 ((22cos(w) + 22cos(2t;) +9cos(3f) + 14cos(4t;) 
+8cos(5 v) - 5cos(6w) + ll)tan3(|) - 27v^)v)/{ cos{v) - 1) 
+{4817v^ + 3120?;2 

+4410) cos(?;) - 2205cos(2 v) - 2205))/((cos(v) - 1)^) 
blnum = -1350^4 csciO(^) + 2251;^ {59v^ + 12) csc^^) 
-5(7361 + 3696^2 - 540) csc^d) + 6(4817 + 3120^^ 
-7710) csc^d) + 264780 CSC 2(|) + 30720(3t;2 - 7)cos(i;) 
+15360 {v^ - 2)cos(2t>) +30fcot(|) ( ^ (-2436 cos (u) 
+1921 cos(2?;) + 1595) csc6( |) + 8322) + 10(5632?;^ 
-1024(8 cos(i') + 23) sm{v)v 
^ 6(65cos(i;) + 66)tan(|)i; _ 36 _ ^ 

cos(t;) + 1 cos(t;) + 1 

blnum = 12150?;^ csciO(§) - 675^2 (191t;2 + 36) csc^d) 

+15 (28499?;'' + 12432i;2 - 1620) csc^d) - 36450i; cot(|) csc^d) 

-12(44896^'^ + 26400i;2 - 35955) esc ^(f) 

+381960?; cot(|) csc'^d) + 12(19268?;-^ + 12480?;^ 

-219885) csc^d) - 817950?; cot(f ) csc^d) - 540 sec^^) 

-5120(73?;2 - 831) - 61440(9t;2 - 40) cos (v) 

-61440(3t;2 - 8) cos(2?;) - 10240(?;2 - 3) cos(3t>) 

-1512900?; cot(f ) + 1884160?; sin (?;) + 532480?; sin(2?;) 

+30720?; sin(3?;) + 90?; sec 2(|)tan(|) + 7740?; tan (|) 

blnum = -4050?;4 csciO(l) + 675?;2 (67?;^ + 12) csc^d ) 

-15(11113 + 4464?;2 - 540) csc6(|) + 12150?; cot (|) csc^^) 

+6 (43259?;^ + 24000?;^ - 24570) csc^(§) - 133920?; cot(|) csc^d 

+347250?; cot(f ) esc 2(|) + 8(4817f'^ + 15840?;2 - 197280) 

+7680(19?;2 - 125) cos(?;) + 30720(2?>2 - 7) cos(2?;) 

+7680(?;2 - 3) cos(3?;) + 292500?; cot (|) - 555520?; sm(?;) 

-194560?; sin(2?;) - 23040 ?;sin(3?;) - 30?; sec^d) tan (f) 

0,0.^ /■"x 72(4817?;'' + 3120?;^ - 26165) 360 
-3180?;tan(-) H — ^ + 



'2 cos (?;) — 1 cos(?') + 1 

H,num = 28350?;^ csc^^d) - 14175?;2 {23v^ + 4) csc^d) 
+45(28269 + 10864?;^ - 1260) csc%) - 85050?; cot(|) csc^(|) 
-12 (182576?;^^ + 98160?;^ - 87255) csc^(f ) 
+965160?; cot(|) csc^(| ) + 180(9634?;^^ + 6240?;^ 
-38111) CSC 2(|) - 2753430?; cot (f) csc^^) - 256(2011?;^ 
+3870?;2 - 45810) - 15360(61?;2 - 479) cos(?;) - 30720(13^2 
-56) cos(2?;) - 76800(?;2 - 3) cos(3?;) - 1054260?; cot (|) 

+3537920?; sin(?;) + 1372160?' sm(2 v) + 230400?; sin(3?;) 

o ,v , ,v , ,v , 360 

-30?;sec2(-)tan (-) - 2100?;tan(-) H — 

^2^ ^2^ ^2^ cos(?;) + l 



QT12_JMC.tex; 15/11/2008; 10:30; p 



High Order Multistep Methods with Improved Phase-Lag Characteristics 



23 



b\num = 930481;'' + 156160f2 + 1042425 csc^d) 
+i5360(9?;2 - 74)cos(i;) - 268800 cos(2v) - 38400 cos(3i;) 
-(9(7560 csc^ (I) + 14963)w^ + 1878401;^ - 4(237663?;'' + 77080t;2 
+298620) cos {v) + 5((34095t;^ - 4000?;2 + 63126) cos(2t;) 
-2v {v{UAhlv'^ - 5200) cos(3t;) + 8(5t;(140 cos(4 v) - 16cos(5i;) 
-lcos{<ov) + 2cos(7t;)) - 666sin(v) + 600sin(2 v) - 298sin(3i;) 
+238sm(4D) + 40sm(5t;) + 80sin(6t;) - 30 sin {7v)) 
+36(sec2(|) + 94)tan( §))) 



270 

+878850) /((cos (v) - if) - 1793280 

cos(t;) + 1 

bldenom = 5120^;^ 

bldenom = 2560i;4 

bldenom = 5120^;^ 
bldenom = 640t;^ 

bldenom = 2560^;'' 

bldenom = 320t;^ 



(37) 



2 _ 90987349 _ 16301796103 v'^ 1568734969 v"^ 

^'^ ~ 53222400 ~ 96864768000 ^ 232475443200 
112423833619 v'^ 128168340031 

"889218570240000 ~ 189225711747072000 " " ' 
2 _ 114798419 16301796103 _ 10540703911 

'^'^ ~ ~ 26611200 ^ 9686476800 „ 44706816000 
871935134531 463103326062011 

^ 44460928512000 ~ 473064279367680000 ^ ' ' ' 

^2 _ 270875723 16301796103 639312597971 v"^ 

'^'2 ~ 17740800 ~ 2152550400 ^ 387459072000 
2472777112313 58110459403753 

11856247603200 ^ 3185618042880000 ~ ' ' ' 
,P 67855831 16301796103 89147839889 
hi. = - + 



^'^ 2217600 , 807206400 16144128000 
1282661397349 3541527647083319 

1482030950400 39422023280640000 " ' 
2 _ 50277247 _ 16301796103 300047111873 

'^'^ ~ 985600 ~ ^ 461260800 27675648000 
18727460555953 4749826752625121 

9880206336000 ^ 22526870446080000 
^2 _ _ 253491379 16301796103 _ 37309797113 

^'^ ~ 4435200 384384000 2767564800 

90333884682737 _ 69221174089601 

^ 37050773760000 250298560512000 



(38) 
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Method PF-D3: 

blnum = -108?;5 csciO(l) + 9v^ {5Sv^ + 12) csc^d ) 

+486f CSC 6(|) - 2484i; csc^(§) - 6U'iv{Av'^ - 9) cos(z;) 
+3cot(|)(99i;2 csc^d) + 9(16 - 9v^) csc^d) - 8(113i;2 
+276) CSC 2(1) _ 11266'y2 + 10032) - 6144(?;=* -6v + 4sin(u)) 
+2^2(26624 sin (t;) - 5635 tan (§)) + 5520tan(|) 

-9 t;scc4(|)(?;tan(|) -6) -6sec2(|)((61 i;^ + 24) tan(|) 

32508 V 

-363v) + -— — 

cosfw) — 1 

blnum = csciO(l) - 451;^ (71v^ + 12) csc^d) 

+i8i;(212 v"^ + 48f2 - 135) csc^Jl) - 14851;^ cot (|) csc^d) 
+159847; CSC ^(f ) + 27(131^2 - 80) cot (|) csc^(|) 
+55782t; esc 2(|) + 6(17777;2 ^ 5952) cot (|) csc^(f ) 
+54z;sec^ (|) + 6144?;(5t;2 - 39) + 6144f (8 - 33) cos(i;) 
+18432t;(t;2 - 4) cos(2i;) + 6 (20035^;^ - 31704) cot(|) 
-20480 (7^2 - 6) sin(u) + 3072(12 - 19^^) sm(2 v) 
-9i;2sec4(|)tan(|) - 12(23i;2 + 12) sec^d) tan (f) 

+14(456 -337.2) tan (^) + ^^ 

2 cos(t;J + 1 

bl^^„^ = -1620?/'^csciO(|) + 135?;^ (85.^ + 12)cscS(f) 
-6.(3872 v"^ + 768.2 _ ^215) csc*^(|) + 4455.2 cot (§) csc^d) 
+36. (424."^ + 96.2 _ ^^ggg^ csc'^d) + 9 (720 
-1781.2) cot(f ) csc''( f ) - 96642. csc2d) - 24 (685.2 
+4716) cot(|) csc2( I) - 54. sec^(|) + 1024. (606 - 67.2) 
-49152.(2.2 - 13) cos(.) - 49152 .(.2 - 5) cos(2.) 
+2048.(21 - 4.2) cos(3 .) + 6(111720 - 46199.2) cot(|) 
+32768(11.2 - 15) sin(.) + 4096(41.2 - 36) sin(2.) 
+4096(7.2 - 6) sin(3.) + 9.2 sec 4(§) tan(|) + 6(61 .2 

+24) sec2(|) tan(|) + 10(1055.2 - 552) tan(|) - 

bl^^^ = 1620.5 csc^°(|) - 135.3 (95.2 ^ ^^2) csc^(f ) 

+18. (1844.'^ + 336.2 _ 4Q5) csc^jf) - 4455.2 cot(|) csc^d ) 

-216.(169.^ + 36.2 _ 288) esc ^{j) + 27(737.2 

-240) cot (f ) csc^d) + 18. (848.^+ 192.2 + 2229) csc2(|) 

+18 (35.2 + 6528) cot (f) CSC ^ (^^) - 54^ sec^ {^) + 384v (133.2 

-1350) + 768.(112.2 _ 843^ ^_ 3072^7(13.2 - 81) cos(2.) 

+1536.(8.2 - 45) cos(3.) + 768.(.2 - 6) cos(4.) + 54(4447 .2 

-13800) cot(f ) + 512(1146 - 659.2 ^ gj^^^^) _ 512(289.2 

-348) sin(2.) + 1024 (42 - 43.2) sin(3.) + 256(12 

, , , o A,v. 9 552 .2 + 288 

-11.2 4 ^ ^ 9^2 gec4 _ + 52587;2 + _L 

^ ^ ^ ^ ^2^ cos(.) + l 

^, , 3276. 

-6384) tan(-) - 



^2' cos(.) + l 



(39) 
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blnum = -11340i;5 csciO(l) + 9A5v^ {lOlv^ + 12) csc^d) 
-18v (15472^4 + 2688?;2 - 2835) csc^d) 

+31185f 2 cot (|) CSC ^(f ) + 108t; (35001;'^ + 720^;^ - 4263) csc^(f ) 

+189(240 - 823t;2) cot (|) csc^^) + 864(87 

-973) cot(l) csc^d ) + 54?; sec^(^) + 384?;(159?;4 - 796 

+8568) + 3072^(1395 - 172?;^) cos(u) - 24576v (11?;^ 

-75) cos(27;) + 3072t>(171 - 28f 2) cos {3v) 

-12288i;(i;2 - 6) cos(4i;) + 18(307848 - 88943?;^ ) cot(|) 

+2048(1055?;2 - 2154) sm{v) + 2048(505?;^ - 708) sin (2?;) 

+2048(155 v"^ - 174)sin(3u) +4096(llz;2 - 12)sin(4i') 

-9 t;2 sec^d) tan(|) - 6 {61v^ + 24) sec2(|) tan (|) + 10(552 

2x 36?;(256(537;2 + 12) u^+ 1479) 4356t; 
-1019?;^) tan (-) + ^ ——^ ^ + — 

^ ^2' cos(t>)-l cos(i;) + l 

blnum = 2268v^ csc^°(f ) - 189^3 (103^2 _^ 12) csc*(f ) 

+hv (9820t;^ + 1680t;2 - 1701) csc^d) - 6237i;2 cot(|) csc^d ) 

-144u(590i;'^ + 1201;^ - 651) esc + 63(511?;^ 

-144) cot (f ) csc''(|) + 18t; (3392t;^+ 768v^ - 677) csc^d) 

+6 (28224 - 3463v^) cot(f ) esc ^(f) + 54t; sec^d) 

-128u {Ulv'^ - 460t;2 + 5016) + 3072i;(32i;2 - 269) cos (v) 

+6144'y(9f;2 - 62) cos(2d) + 2048?;(8?;2 - 51) cos(3i;) 

+3072i;(t;2 - 6) cos(4w) + 1002(317 v'^ - 1128) cot(f ) 

-6144(67f 2 - 146) sin(t>) - 3072(69t;2 - 100) sin(2t;) 

-12288(5 - 6) sin(3f ) + 1024(12 - llv^) sin(4 v) 

_(9^2 gg^4( 543g^2 552 i;^ + 288 _ v 

^ ^2^ cos(i;) + l ^ ^2^ 

3276^ 

cos(f ) + 1 
bidenom = 3072^^ 

bldenom = 1536?;5 

bUnom = 1024?;5 
bldenom = ^^iv' 
bUenom = 
bldenom = 256?^^ 



(40) 
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3 _ 90987349 _ 16301796103 1273143229t;^ 

^'^ ~ 53222400 ~ 72648576000 ^ 193729536000 
421960559 124506164597657 



1221454080000 5676771352412160000 
3 _ 114798419 16301796103 _ 38969308351^^ 

"^'^ ~ ~ 26611200 ^ 7264857600 , 96864768000 
5604849953 'if _ 5242236508523657 

^ 158789030400 ~ 2838385676206080000 

^3 _ 270875723 16301796103 578948920601 

~ 17740800 , 1614412800 ^ 193729536000 
733132845253 3064860132866873 ^ 

1482030950400 ^ 57341124771840000 
3 __ 67855831 16301796103 _ 247322293877 

^'^ ~ ~ 2217600 605404800 24216192000 

844539338551 80127849554001773 

^ 370507737600 236532139683840000 ^ "' 
3 _ 50277247 _ 16301796103 279925886083 

'^'^ ~ 985600 ~ 345945600 ^ 13837824000 
644885438797 122198747647321067?;^ 

123502579200 ^ 135161222676480000 , 
^3 __ 253491379 16301796103 1220113637^^ 

^'^ ~ 4435200 ^ 288288000 „ 48384000 
55089606671 1634963205351829 



(41) 



8096760000 1325110026240000 



+ ... 
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Method PF-D4 : 

6f „„„ = 216v^ csc^^d) - SlOv'^ CSC ^(f ) - 54(157^2 
+80) CSC 4(f) + 36(1840 - 2257i;2) esc 2(|) - 12288(20?;'^ 
-71^2 + 20) cos (t;) - 3t;cot(§)((4682t>2 + 9csc 2(|)(5csc2(§) v'^ 
+5Av^ + 96) + 5856) csc^d) + 2520 {Slv"^ - 48)) 
-45t;2 sec^d ) (vtan(|) - 6) + 56v(1024 {llv'^ - 12) sm{v) 
+135(48 -31f 2) tan (|)) - 18sec4(|) (-373t;2 + 3(21i;2 
+16)tan(f) t;-80) - 2 sec2(|)(-39246 v'^ 
+{6733v'^ + 6384) tan(|)z; + 28320) 

binum = -W80v^ csci°(f ) + 2160?;'' csc^d) + 4050^^ csc^d) 

+54(617i;2 + 400) csc^d) + 293508i;2 csc^d) + 12288(20 -y^ 

-71^2 + 20) cos(u) + 24576(10t;4 - 59t;2 + 20) cos(2t;) 

+3V cot ( f ) ( ( 16978?;2 + 9 csc^ ( f ) (25 csc^ ( § ) z;2 + 2 14^^ + 480) 

+20064) csc^d) + 120(2341t;2 - 5136)) + (72(4096 

-31783z;2 - 80(257cos(w) + 251) csc2(t>) + 1024(4t>4 - 35 v"^ 

+20)cos(t;) +25600))/(cos(v) + 1) - v{135v {vtan(^) 

-6) sec 6(|) + 18((161u2 + 144) tan (f) - 951i;) sec^(f ) 

+18(1407?;^ + 1616)tan(|) sec^d) + 8(7168(lli;2 _ 12) sm(t;) 

+1024(107^2 - 156)sin(2?;) + 15(3217 v'^ - 9168) tan(f ))) 

binum = 9720?;^ csc^O(|) - 34560^^ csc^d) + 270z;2(128i;4 

-135) CSC 6(1 ) - 6075t;''^ cot(|) csc'^d) - 162(1459?;^ 

+1200) CSC 4(1) - 54t;(767f 2 + 2160) cot( f ) csc^d) 

+108 (32720 - 18831i;2)csc2(|) - 234?; (1511?;2 

+1488) cot(f ) CSC 2(|) + 3510?;2 sec^d ) + 90(701?;2 

+208) sec^d) - 12288(98f4 - 825?;2 + 540) - 12288(220?;^ 

-1717 + 940) cos(?;) - 98304(10?;4 - 59i;2 + 20) cos(2 v) 

-24576(20?;4 - 147?;2 + 60) cos(3i;) + 360t; (45744 

-16567i;2) cotd) + 8192z; (1279t;2 - 2868) sin(?;) + 32768t;(107i;2 

-156) sin(2i;) + 147456f (13f 2 - 24) sin(3t;) - 7;(585?;2 scc^d) 

+(42991t;2 + (32281i;2 + 46128 ) cos(t;) + 57360) sec^d) 

/ 2 ^^ ^v, 24(30457^2 -69040) 
+120 (8741t>2 - 38928)) tan(-) + ^.q^^(^) ^ ^ 

bi,num = -3240?;^ csc^Od) + 15120?;^ csc^d) - 810z;2(32w4 
-15) CSC 6(|) + 2025?;3cot(J) csc^d) + 54(320t;<' + 1179i;2 
+1200) csc^d) + 54v(209v^+ 720) cot (f) csc^d) 
+108 (5237t;2 - 11440) csc2d) + 18^; (5461?;2 
+3888) cot(|) esc 2(1) + 270t;2 sec^d) + 90 (41i;2 + 16) sec^d) 
+76800(5 v"^- 54^2 + 48) + 6144(100^;^ - 819?;2 + 480) cos (v) 
+73728(5z;4 - 39z;2 + 20) cos(2?;) + 6144(20 v"^ - 147z;2 
+60) cos(3t;) + 6144(5?;4 - 4l7;2 _^ 20) cos(47;) + 360?;(4705t;2 
-14928) cotd) - 36864i;(67i;2 - 154) sin(t>) - 86016w(17 v'^ 
-36) sin(2t;) - 36864?;(13i;2 - 24) sin(3 v) - 2048t;(61i;2 

90 

-132) sin(4t;) - v (((1087 + — + (457 

cos(f ) + 1 

+1776) cos(i;) + 2640) sec4(^) + 120(349t;2 

24(349^2 - 5680) 
-2832))tan(-) + ^ ^ 



'2^ cos(u) + l 



(42) 
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bl^,^^^„ = 22680i;6 csc^°(|) - 120960i;6 csc^d) + 4050t;2 (64^;^ 
-21) CSC - 14175^3 cot(§) csc^d) ~ 54(5120?;^ + 7077t;2 
+8400) csc''(^) - 378i;(181i;2 + 720) cot(|) csc'^d) + 108 (1280t; 
-33131?;2 + 82320) csc^^) - 378?;(1643v2 + 784) cot(|) csc^d) 
-1890f2 sec'^d) - 126(277t;2 + 80) sec^d ) - 15360(160t;'^ 
-1833f 2 + 1740) - 30720(120 - 1079t;2 + 740) cos(t;) 
-6144(400?;'^ - 3271^^ + 1780) cos(2u) - 6144(140?;^ - 1113^^ 
+540) cos (3 v) - 49152(5t;^ - Alv'^ + 20) 008(41;) - 6144(4^;^ 
-35 iP + 20) cos (5?;) + 360t;( 102576 - 30047f 2) cot (§) 
+20480i;(755i;2 - 1956) sin (v) + 4096i;(2429i;2 - 5412) sm(2t;) 
+12288?;(281 v'^ - 588) sin(3?;) + 16384i;(61i;2 - 132) sm(4 v) 
+20480^(5^2 - 12)sin(5v) + v (((25429 + 33^^)^^ + (19507 v 

+26256)008(1;) + 32304) seo'^(|) + 120(53471;^ - 20976)) tan(|) 
885120 - 443976i;2 

0O8(v) + 1 

binum = -13608t;6 csc^^{^) + 756001;^ cso^d) + 270^;^ (189 

-MOv^) 080 6(|) + 8505^;^ cot(|) oso^d ) + 162(1280i;6 

+1337i;2 + 1680) csc^d) + 378t;(103?;2 + 432) cot(f ) csc^d ) 

+126i; (2863t;2 + 1104) cot(|) esc 2(|) _ 18901;^ sec^d) 

-126(253i;2 + 80) sec^(§) + 3072(181;^ + 400?;^ - 4575^2 + 4500) 

+30720(80 - 753^2 + 540) oo8(u) + 30720(40w^ - 327v^ 

+180) cos(2v) + 6144(100i;^ - 8191;^ + 420) cos(3i;) + 24576 (5u' 

-41^2 + 20) cos(4u) + 6144(4?;"^ - 35^^ + 20) cos(5i;) 

+360^(17495^2 - 61104) cot (f) - 552960^(19^2 - 52) sin (v) 

-552960?;(9w2 - 20)sin(2i;) - 12288i;(203 v'^ - 444)sin(3i;) 

-8192z;(61f2 - 132)sin(4 v) - 20480t;(5f2 - 12)sin(5i;) 

+z;(63(48 (v'^ + 1) + (43?;2 + 48) cos(t;)) sec ^d) + 840(659?;2 

, 57961 + 86592, 

-2928 + )tan - 

^ cos(d) + 1 ^ ^2' 

_^216 (1280i;^ - 19257^2 + 49840) _^ 908160 - 325608?;2 

cos {v) — 1 cos(t;) + 1 

b^nom = 24576.;6 

bUnom = 12288t;6 

bUnom = 24576t;6 

bXaenom = 3072^;^ 

bUnom = 12288.;6 

bUnom = 6144^;^ 



,6 



(43) 
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4 _ 90987349 _ 16301796103 v"^ 451826777 v"^ 

^'^ ~ 53222400 ~ 58118860800 ^ 87178291200 
567827928197 7456559915267 



711374856192000 87334943883264000 
4 _ 114798419 16301796103 _ 53423656079 

~ ~ 26611200 ^ 5811886080 , 87178291200 
3768047264957 1807362562020419^;^ 

71137485619200 ~ 567677135241216000 ^ ' ' ' 

^4 _ 270875723 16301796103 137191770479 

^'^ ~ 17740800 ~ 1291530240 ^ 29059430400 
6482264498201 177746918811101 



6774998630400 1557413265408000 , (aa^ 

4 __ 67855831 16301796103 _ 4393734833?;^ ^ ' 

^'"^ ~ ~ 2217600 't 484323840 269068800 

3998343544387 42449166474915191 

^ 846874828800 47306427936768000 " ' 

4 _ 50277247 _ 16301796103 9638045471?;'' 

'^'^ ~ 985600 ~ , 276756480 296524800 
87748735805873 70694795981650949 _ 

7904165068800 ^ 27032244535296000 , 
^4 _ 253491379 16301796103 84219941177 

"^'^ ~ 4435200 \ 230630400 2075673600 
863717128805537 _ 1830923526668191 

^ 59281238016000 500597121024000 ^ ' ' ' 
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Method PF-D5: 



bi,num = - cos^(l) csc^^{v) siii ^ (|) (60(-5i;^ + 3v^ + 3sm{v) 
—6 siii(3?;) + 3sin(4u) + 3siii(5f) — 6sin(6f) + 3sm(8f)) 
+i?;(15(200z;'^ + 283v'^ - 100) cos (v) + 300(2 cos(2t;) + 6 cos(3t;) 
— 10 008(4?;) + cos(5f) + 12 cos {6v) — 2cos(7f) — 4cos(8t;)) 
+t;(t;(15(366 cos(3 v) + 686 cos(4t;) - 503 cos(5t;) - 480 cos(6t;) 
+142cos(7w) + 116 cos(8i;)) + 2(15(72?;2 - 155) cos(2w) 
+i;(-60 i;(45cos(3i;) + 24cos(47;) - 25cos(5z;) - 14cos(6i;) 
+5cos(7 v) + 3cos(8t;)) + 405cot(|) + 357sin(i;) - 2530 sm{2v) 
+4516sin(3t;) + 3752sin(4z;) - 3373sin(5?;) - 2324 sin {6v) 
+770sm(7f) + 522sin(8w) +315tan(f) ))) - 60(49 sin(t>) 
-44sin(2t;) - 6sm(3i;) + 131 sin(4i;) - 67sin(5i;) - 114sin(6i;) 
+28 sin(7v) + 31 sin(8f ))))) 

blnum = csc^d) sec9(|)(4(3875?;^ - 519?;^ - 90) cos(u) 

-36(147 - 76v^ - 10cos(2u) + 20cos(4z;) - 20cos(5t;) 

+10cos(7i;) - 20 cos(8t;) + 10 cos(10z;)) + t;(t;(2(5117 v'^ 

-4254) cos (2z;) + 12(516 cos (3?;) + 482 cos(4t;) - 158 cos {5v) 

+348cos(6?;) - 317cos(7t;) - 486cos(8i;) + 132cos(9?;) 

+137 cos(10t;)) +t;(3(-12sm(t>) - 3789sm(2t;) +3900sin(3 v) 

+2058sin(4t;) + 708sm(5u) + 1632sin(6i;) - 2142sin(7 v) 

-1792 sin(8i;) + 614 sin(9u) +461 sm(10t;)) 

-2i;(-393216 u(15cos(2i;) + 17) sin'^d) cos ^(f) + 4476 cos (3z;) 

+1918cos(4t;) + 2288 cos {5v) + 1242cos(6w) - 2476cos(7i;) 

-1400 cos(8^;) + 582 cos(9 v) + 351 cos(10?;)))) 

-15360COS 3(f )(186cos(w) + 156cos(2t;) + 132cos(3 v) 

+78cos(4?;) +48cos(5t;) + 19cos(6t;) + 98) sin 5(|))) 

blnum = -cos6(|)csc^4(^;)sin 8(|)(7290cot(§) 

+5670tan(f )7;'^ + 32(56^^ - 93) sin(9f )7;2 + 12(250v^ - 751v'^ 

+WO)v + 3(2140 + 13737?;2 - 2220) cos(D)t> - 6(260f'' 

-591^2 + 540) cos(2t;)t; + 6(-1702«'i + 2971z;2 + 1100) cos(3t;) v 

-6(584?;'' - 3933v^ + 1020) cos(4t;)t; + 15(236 v'^ - I395v^ 

+228) cos(5?;)t; + 12(182^;^ - 1121^;^ + 800) cos(6v)t; + 12(61?;'' 

-201f 2 - 180) cos(7t;)t; + 12 (lOv^ - 171t;2 + 60) cos(8t')i; 

-60(8?;^ - 53v^ + 20) cos(9f )f - 24(10?;^ - 93?;^ + 90) cos(10?;)i; 

+6 (1327?;^^ - 974z;2 - 60) sin(t;) - 4(1295?;^ + 1098 v"^ 

-180) sin(2i;) + 8(3403;;'i + 579^2 - 270) sin(3 v) + 4(3256?/^ 

-4557t;2 + 90) sin(4t;) + 2(-6097 + 8058t;2 + 180) sin(5i;) 

-8(935?;'' - 1929^2 + 270) sin(6?;) - 4(593^" + UGv"^ 

-180) sin(7i;) - 4 (u^ - 3)(199?;2 + 30) sin(8?;) + 8(127 - 363?;^ 

+90) sin(10'i;)) 

(45) 
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binum = 5670i;3 sec«(|) - 945^'' tan (f ) sec^(|) + 270v {257v^ 

+160) sec6(|) - 1080^;2 (ii„2 ^ ^^g^ tan(f ) sec + 4320t;(51?;2 
-46) sec ^(|) - 18(2267i;4 + 9816?;2 + 2880) tan(|) sec^(|) 
+4 (-eieSlv'i + 89664^2 + 570240) tan(f ) sec2(^) 
+7290i;3csc 6(|) + I215w^ cot(|) csc^^) + 78732^;^ csc^(|) 
+81t;2(l67t;2 + 240) cot(|) csc^d ) + 756t;(1127t;2 
-2640) CSC 2(|) + 9(16297?;^ + 2592^2 - 17280) cot(|) csc^d) 
+46080?; (30?;^ - 443?;2 + 980) + 73728?;(101?;2 - 370) cos {v) 
+36864z;(40z;^ - 431?;^ + 580) cos(2u) + 73728^; (13t;2 
-40) cos(37;) + 6144t;(30t;^ - 307t;2 + 340) cos(47;) + 9(342973?;'' 
-1171680^2 + 328320) cot (|) + 73728(17?;i - 271?;^ 
+210) sin {v) - 98304(67?;'' - 243?;^ + 90) sin(2t;) 
+16384(10 - 147?; 2 + 90) sin(3?') - 24576(33?;^ - 107?;2 
+30) sin(4?;) + (-3937883?;'' + 28223520?;2 - 21006720) tan (|) 
24?;(115673?;2 - 893360) 
cos(?;) + 1 

blnum = 39690?;''^ scc^d ) - 6615?;'' tan (|) sec^d) + 1890?- (239?;2 
+160) sec6(|) - 1890?;2 (41?;2 + 72) tan(|) sec 6(|) + 756?'(1471?;2 
-2160) sec ^(|) - 378(561?;^ + 2968?;2 + 960) tan (|) sec^d ) 
+84?; (99661?'2 - 859120) sec2(|) - 56 (26371?;^ - 61404?;2 
-289440) tan(f) scc2(|) - 51030?'3 csc^d) 
-8505?;^ cot(|) csc^d ) - 503496?;^ csc''(|) - 1701?;2(51 ?;2 
+80)cot(f)csc^(f) - 756?;(7131?;2 - 18320) csc2(|) 
-189(4913?;^ + 256?;2 - 5760) cot(|) csc2(|) - 92160?;(373?;2 
-2140) - 245760?; (60?;" - 767?;2 + 1400) cos(?;) - 12288?; (1829?;2 
-6620) cos(2?;) - 884736?; (5?;" - 54?;2 + 70) cos(3?;) - 24576?;(61?;2 
-220) cos(4?)) - 49152?; (6?;" - 65?;2 + 80) cos(5?') - 9(2193517?;" 
-8002080 ?;2 + 2378880) cot(§) + 245760(287 ?>'' - 1321?;2 
+690) sin(?>) - 24576(151?;" - 2446?;2 + 1860) sin(2?;) 
+16384(1211?;" - 4317?;2 + 1530)sin(3 ?;) - 49152(5?;" - 82?;2 
+60) sm(4?^) + 49152(27 ?'" - 95?'2 _^ 39) sin(5?;) + (-24636373?;" 
+185682720 ?;2 - 157806720) tan(|) 
blnum = 119070?;=^ sec^d) - 19845?;" tan (|)sec8(|) 
+5670?; (233W2 + 160) sec'^d) - 45360^2 (5,^2 _^ tan(|) sec ^(|) 
+30240?- (97?'2 - 170) sec "(f) - 126(4501?;" + 25800?;2 
+8640) tan(|) sec"(f ) + 28 (-151573?;" + 398400?;2 
+1745280) tan(f ) sec2(|) + 153090?;^ csc^d) 
+25515?-" cot(f ) csc^d) + 1462860?;3 csc"(f ) + 2835?;2 (89?;2 
+144)cot(§) csc"(f) + 3780?;(4139?;2 - 10960) csc2(|) 
+63(42779?;" + 480?;2 - 51840) cot(|) csc2(|) + 6 1440?- (400?;" 
-6273?;2 + 15540) + 184320?;(737 ?;2 - 3260) cos(?;) 
+184320?;(150?-" - 1739?;2 + 2660) cos(2?;) + 368640?; (61?;2 
-220) cos(3?;) + 61440?; (72?;" - 779?;2 + 980) cos(4?;) 
+614400?;(?;2 - 4) cos(5?;) + 30720?;(4?;" - 45?;2 + 60) cos(6?;) 
+27 (2129819?;" - 7933600?;2 + 2405760) cot(|) + 122880(173?;" 
-3222?;2 + 3060) sm(?;) - 491520 (259?;" - 1035?'2 + 450) sm(2?;) 
+8192(451?;" - 7365 ?;2 + 5490) sin(3?;) - 24576(809?;" - 2860?;2 

(46) 



QTi2_JMC.tex; 15/11/2008; 10:30; p. 31 



32 



D.S. Vlachos, Z.A. Anastassi, T.E. Simos 



+960) sm{Av) + 49152(2?;'' - 35v'^ + 30) sm{5v) - 4096 (137?;'' 
-510v^ + 180) sin(6?;) + (-71232167 v"^ + 545068320t;2 

.o.^.o^o N .V. 4200i;(11495t;2 - 101648) 

-484318080 tan(-) + ^ 

^ ^2^ cos(?;) + l 

^'I^eno^ = 393216.^ (47) 

^3,denom = 

^!,deno™ = 6144t;7 
&i,denon^ = 24576?;^ 
bUnom = 61440?;^ 



, . 90987349 16301796103 3000167803 



53222400 , 48432384000 1162377216000 
2795461105681 423149542577683 



1778437140480000 1892257117470720000 

5 _ 114798419 16301796103 9164631311 v'^ 

'^'2 ~ ~ 26611200 4843238400 10567065600 ^ 
1149042063431 _ 1015140295561079i;^ _ 

16167610368000 ~ 189225711747072000 ~ ' ' ' 

5 _ 270875723 16301796103 106131595741 v"^ 

^'^ ~ 17740800 ~ 1076275200 ^ 15498362880 
193910381341843 12701953007551 

118562476032000 ^ 62296530616320 ' ' ' , (Aii) 
^5 _ 67855831 16301796103 _ 15415020883 v"^ ^ ^ 

^'^ ~ 2217600 t 403603200 645765120 ^ 

125531790304181 361081038637727 



+ ■ 



14820309504000 185515403673600 
5 _ 50277247 16301796103 263828905451 v'^ 

^'^ ~ 985600 ~ „ 230630400 5535129600 
400330948542829 54441481708734389 

+ 



19760412672000 9010748178432000 
5 _ _ 253491379 16301796103 _ 117727186937 v"^ 

^'^ ^ 4435200 t 192192000 1976832000 

3958504514434801 7164842887862063 

+ ... 



148203095040000 834328535040000 
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5— aiiis 

Figure 1. The stability region (s-v plane) of the classical Quinlan-Tremaine 12-step 
method 
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Figure 2. The stability region of the methods PF-D0,PF-D1, PF-D2, PF-D3, PF-D4 
and PF-D5 (from left to right and from top to bottom) 



-■-MJJTI 2 

-»-M_QT1 2_PF 
^- M_QT1 2_PF_D1 
Tt-M_QT1 2_PF_D2 
-•-M_QT1 2_PF_D3 
< M_QT1 2_PF_D4 
-H-M QT12 PF D5 




Step Size 

Figure 3. The accuracy (digits) of the new methods compared to the classical one 
for the Schrodinger equation {E = 989.701916) 
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step Size 



-♦-M_QT12_EF 

^^M_QT12_EF_D1 
-*-M_QT12_EF_D2 
-►-M_QT12_EF_D3 
'«-M_QT12_EF_D4 
-M-M QT12 EF 06 




Figure 4- The accuracy (digits) of the new methods compared to the classical one 
for the Schrodinger equation (E = 341.495874) 



M_QT12 

M_QT12_EF 

M_QT12_EF_D1 

M_QT12_EF_D2 

M_QT12_EF_D3 

M_QT12_EF_D4 

M QT12 EF 05 




Step Size 



Figure 5. The accuracy (digits) of the new methods compared to the classical one 
for the Schrodinger equation {E — 163.215341) 
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